* Load dataset
use "Yourpath\080224_combineddata_clean.dta", clear

* Save backup file to prevent data loss
save "Yourpath\100424_unconstrainedIDL.dta", replace

* Keep observations where phase == 2
keep if phase == 2

* Standardize variables
egen zmobility = std(mobility)
egen zlgwdp = std(lgwdp)

* Sort and generate means
sort fips_code time
by fips_code: egen mzmobility = mean(zmobility)
by fips_code: egen mzlgwdp = mean(zlgwdp)
gen dzmobility = zmobility - mzmobility
gen dzlgwdp = zlgwdp - mzlgwdp
tab cyclic, gen(dcyclic)

* Generate 'used' variable
gen used = 1 if time >= 58
replace used = 0 if time < 58

sort fips_code used time

* Generate lagged variables for zlgwdp (1 to 20 lags)
foreach i of numlist 1/20 {
    by fips_code: gen lag`i'_zlgwdp = zlgwdp[_n - `i']
}

* Generate lagged variables for zmobility (1 to 4 lags)
foreach i of numlist 1/4 {
    by fips_code: gen lag`i'_zmobility = zmobility[_n - `i']
}

* Generate mean-adjusted zmobility
by fips_code used: egen mzmobility2 = mean(zmobility)
gen dzmobility2 = zmobility - mzmobility2

* Get unique list of states
levelsof mst, local(states)

* Define output Excel file path
local output_file = "Yourpath\100424_unconstrained_IDLcoefficient.xlsx"

* Initialize Excel file and specify Summary and Coefficients sheets
putexcel set "`output_file'", replace

* ---------------------------------------------------------------------------
* Sheet for Summary
* ---------------------------------------------------------------------------
putexcel set "`output_file'", modify sheet("Summary")
putexcel A1 = "State" B1 = "Model Type" C1 = "Overall Impact" D1 = "R-squared"

* Initialize row number for the "Summary" sheet
local row = 2

* ---------------------------------------------------------------------------
* Sheet for Coefficients
* ---------------------------------------------------------------------------
putexcel set "`output_file'", modify sheet("Coefficients")
putexcel A1 = "State" B1 = "Lag" C1 = "Coefficient"

* Initialize row number for the "Coefficients" sheet
local coef_row = 2

* ---------------------------------------------------------------------------
* Run unconstrained model at the National level
* ---------------------------------------------------------------------------
di "Running Unconstrained model at the National level"

preserve

* Ensure there are observations at the national level
count
if r(N) == 0 {
    di "No observations at the national level. Exiting..."
    restore
    exit
}

* Step 1: Run unconstrained IDL model
xtreg zmobility lag1_zlgwdp-lag20_zlgwdp spring fall summer i.cyclic, fe coeflegend

* Capture coefficient names
local coefnames : colnames e(b)

* Step 2: Calculate Overall Impact
local overall_impact_unconstrained = 0

* Initialize local macros to store all lag coefficients
foreach i of numlist 1/20 {
    local varname = "lag`i'_zlgwdp"
    local exact_name = ""
    foreach cname of local coefnames {
        if strpos("`cname'", "`varname'") {
            local exact_name = "`cname'"
            break
        }
    }
    if "`exact_name'" != "" {
        local coef_lag`i' = _b["`exact_name'"]
        local overall_impact_unconstrained = `overall_impact_unconstrained' + `coef_lag`i''

        * Step 7: Export each lag coefficient to the Coefficients sheet
        putexcel set "`output_file'", modify sheet("Coefficients")
        putexcel A`coef_row' = "National" B`coef_row' = `i' C`coef_row' = `coef_lag`i''
        local coef_row = `coef_row' + 1
    }
}

* Step 3: Generate predicted variable pred_mob_d (based on code 4 logic)
gen pred_mob_d = 0
foreach i of numlist 1/20 {
    replace pred_mob_d = pred_mob_d + (_b["lag`i'_zlgwdp"] * lag`i'_zlgwdp)
}

* Step 4: Run new xtreg regression and get R-squared
xtreg zmobility pred_mob_d spring fall summer i.cyclic, fe

* Step 5: Get within R-squared
scalar r_squared_unconstrained = e(r2_w)

* Step 6: Export national-level Summary results to the "Summary" sheet
putexcel set "`output_file'", modify sheet("Summary")
putexcel A`row' = "National" B`row' = "Unconstrained" ///
    C`row' = `overall_impact_unconstrained' D`row' = r_squared_unconstrained

* Increment row number to write the next result
local row = `row' + 1

restore

* ---------------------------------------------------------------------------
* Run unconstrained models for each state
* ---------------------------------------------------------------------------
foreach state of local states {
    local clean_state = subinstr("`state'", " ", "_", .)

    di "Running Unconstrained model for state: `state'"

    preserve
    keep if mst == "`state'"

    count
    if r(N) == 0 {
        di "No observations for state: `state'. Skipping..."
        restore
        continue
    }

    * Run unconstrained IDL model
    xtreg zmobility lag1_zlgwdp-lag20_zlgwdp spring fall summer i.cyclic, fe coeflegend

    * Capture coefficient names
    local coefnames : colnames e(b)

    * Calculate Overall Impact
    local overall_impact_unconstrained = 0

    * Export state-level lag coefficients to the Coefficients sheet
    foreach i of numlist 1/20 {
        local varname = "lag`i'_zlgwdp"
        local exact_name = ""
        foreach cname of local coefnames {
            if strpos("`cname'", "`varname'") {
                local exact_name = "`cname'"
                break
            }
        }
        if "`exact_name'" != "" {
            local coef_lag`i' = _b["`exact_name'"]
            local overall_impact_unconstrained = `overall_impact_unconstrained' + `coef_lag`i''

            * Export each lag coefficient to the Coefficients sheet
            putexcel set "`output_file'", modify sheet("Coefficients")
            putexcel A`coef_row' = "`clean_state'" B`coef_row' = `i' C`coef_row' = `coef_lag`i''
            local coef_row = `coef_row' + 1
        }
    }

    * Generate predicted variable pred_mob_d
    gen pred_mob_d = 0
    foreach i of numlist 1/20 {
        replace pred_mob_d = pred_mob_d + (_b["lag`i'_zlgwdp"] * lag`i'_zlgwdp)
    }

    * Run new xtreg regression and get R-squared
    xtreg zmobility pred_mob_d spring fall summer i.cyclic, fe

    * Get within R-squared
    scalar r_squared_unconstrained = e(r2_w)

    * Export state-level Summary results to the "Summary" sheet
    putexcel set "`output_file'", modify sheet("Summary")
    putexcel A`row' = "`clean_state'" B`row' = "Unconstrained" ///
        C`row' = `overall_impact_unconstrained' D`row' = r_squared_unconstrained

    * Increment row number to write the next state's Summary result
    local row = `row' + 1

    restore
}

di "Results successfully saved to `output_file'"
